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Abstract 

We show how to modify the Hnear-time construction algorithm for sufRx arrays based on 
induced sorting (Nong et ah, DCC'09) such that it computes the array of longest common prefixes 
(LCP-array) as welL Practical tests show that this outperforms recent LCP-array construction 
algorithms (Gog and Ohlebusch, ALENEX'll). 

1 Introduction 

The suffix array is an important data structure in text indexing. It is used to solve many tasks 
in string processing, from exact and inexact string matching to more involved tasks such as data 
compression, repeat recognition, and text mining. It is also the basic building block for the more 
complex text index called the suffix tree, either indirectly for index construction, or directly when 
dealing with compressed suffix trees. In all of the above applications (possibly apart from exact 
string matching), the suffix array is accompanied by its sister-array, the array of longest common 
prefixes (LCP-array for short). 

Since their introduction in the early 1990's, much research has been devoted to the fast construc- 
tion of suffix arrays. Although it is in principle possible to derive the suffix array from the suffix 
tree, for which linear-time algorithms had already been discovered earlier [22], for reasons of time 
and space the aim was to construct the suffix array directly, without help of the tree. This long 
line of research (see [20] for a good reference) culminated in three linear-time algorithms [11 U 13 1 IT^ . 
However, these algorithms were notorious for being "linear but not fast" [2], as they were slower 
than other non-linear algorithms that had been discovered before and continued to be discovered 
afterwards. 

This un-satisfactory situation (at least for theoretical practioners or practical theoreticians, who 
want linear-time algorithms to perform faster than super-linear ones) changed substantially when in 
2009 a new linear-time algorithm based on induced sorting was presented [H]. A careful implemen- 
tation of this approach due to Yuta Mori led to one of the fastest known suffix array construction 
algorithms, often outperforming all other linear or super-linear implementations. 

Less emphasis has been put on the efficient construction of the LCP-array. Manber and Myers |15j 
mentioned that it can be constructed along with their method for constructing the suffix array, but 
their algorithm ran in O(nlgn) time and performed rather poor in practice. Kasai et al. [12] 
gave an elegant algorithm for constructing the LCP-array in linear time. A few refinements of 
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this algorithm led to improvements in either space [16j or in running time [TO]. However, these 
algorithms could not compete with the carefully tuned algorithms for suffix arrays. This led to 
the odd situation that the rather difficult task of sorting suffixes could be solved faster than the 
seemingly simpler task of computing longest common prefixes. 

This situation changed only recently when a theoretically slow 0{n?) but practically fast LCP- 
array construction algorithm was presented [7] . Their algorithm exploits properties of the Burrows- 
Wheeler- Transformation (BWT) of the text, which must be computed before. The authors of [7] 
also sketch how their approach yields a linear-time algorithm (for constant alphabets, otherwise it 
takes 0{nlga) time). 

Driven by the success of the fast linear-time algorithm based on induced sorting [18] , we show in 
this paper how it can be adapted such that it also induces the LCP-values (Sect. [3]). This results 
in a new linear-time algorithm for constructing LCP-arrays (for integer alphabets). In Sect. H] we 
show that an ad-hoc implementation of the theoretical ideas leads to a fast practical algorithm that 
outperforms all other previous algorithms. An additional advantage of our algorithm is that it does 
not rely on the BWT, and is hence preferable in situations where the BWT is not already present 
(such as compressed suffix arrays not based on the BWT [T7], for example). 

Before detailing our theoretical and practical contributions, in Sect. [2] we first introduce some 
notations, and then review the induced sorting algorithm for suffix arrays. 

2 Previous Work and Concepts 

2.1 Suffix- and LCP-Arrays 

Let T = ti . . .tn he & text consisting of n characters drawn from an ordered alphabet S of size 
a = |$]|. The substring of T ranging from i to j is denoted by Ti j, for 1 < i < j < n. The substring 
Tj,,„ is called the i'th suffix of T and is denoted by Si. As usual, for convenience we assume that 
T ends in a unique character $ which is not present elsewhere in the text, and that $ < a for all 
a G S. 

The suffix array SA[1, n] of T is a permuation of the integers in [1, n] such that 5'sA[i-i] <iex SsA[i] 
for all 1 < i < n. In other words, SA describes the lexicographic order of the suffixes. 

The array LCP of longest common prefixes is based on the suffix array. It holds the lengths of 
the longest common prefixes of lexicographically adjacent suffixes, in symbols: LCP[i] = max{^ > 
1 TsA[i]..SA[i]+£-i=SA[i-i]..SA[i-i]+f-i} for 1 < i < n, and LCP[0] = 0. 

2.2 Constructing Suffix Arrays by Induced Sorting 

As the basis of our new LCP-array construction algorithm is the induced sorting algorithm for 
constructing suffix arrays |18) . we explain this latter algorithm in the following. Induced sorting 
has a venerable history in suffix sorting, see [HKMlET]. Its basic idea is to sort a certain subset of 
suffixes, either directly or recursively, and then use this result to induce the order of the remaining 
suffixes. In the rest of this section, we follow the presentation of Okanohara and Sadakane |19) . 

Definition 1. For 1 < i < n, suffix Si is said to be S-type if Si <iex Si+i, and L-type otherwise. 
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The last suffix is defined to be S-type. 

The type of each suffix can be determined in hnear time by a right-to-left scan of T: first, Sn is 
declared as S-type. Then, for every i from n — 1 to 1, is classified by the following rule: 

Si is S-type iff either ti < tj+i, or ti = tj+i and S'j+i is S-type. 

We further say that an S-suffix Si is of type S* iS Si-i is of type L. 

In SA, all suffixes starting with the same character c € S form a consecutive interval, called the c- 
bucket henceforth. Oberve that in any c-bucket, the L-suffixes precede the S-suffixes. Consequently, 
we can sub-divide buckets into S-type buckets and L-type buckets. 

Now the induced sorting algorithm can be explained as follows: 

1. Sort the S*-suffixes. This step will be explained in more detail below. 

2. Put the sorted S*-suffixes into their corresponding S-buckets, without changing their order. 

3. Induce the order of the L-suffixes by scanning SA from left to right: for every position i in SA, 
if SsAfij-i is L-type, write SA[i] — 1 to the current head of the L-type c-bucket (c = isA[j]-i)) 
and increase the current head of that bucket by one. Note that this step can only induce "to 
the right" (the current head of the c-bucket is larger than i). 

4. Induce the order of the S-suffixes by scanning SA from right to left: for every position i in SA, 
if 5'sA[i]-i is S-type, write SA[i] — 1 to the current end of the S-type c-bucket (c = isA[j]-i)) 
and decrease the current end of that bucket by one. Note that this step can only induce "to 
the left," and might intermingle S-suffixes with S*-suffixes. 

It remains to explain how the S*-suffixes are sorted (step 1 above). To this end, we define: 

Definition 2. An S*-substring is a substring Ti„j with i ^ j of T such that both Si and Sj are 
S*-type, but no suffix in between i and j is also of type S*. 

Let i?2> • • • > Rn' denote these S*-substrings, and a' be the number of different S*-substrings. 
We assign a name Vi G [l,cr'] to any such Ri, such that Vi < vj if Ri <iex Rj and Vi = vj if 
Ri =iex Rj- We then construct a new text T' = vi...Vn' over the alphabet [1,0"'], and build 
the suffix array SA' of T' by applying the inducing sorting algorithm recursively to T' if a' < n' 
(otherwise there is nothing to sort). The crucial property [TH] to observe here is that the order of 
the suffixes in T' is the same as the order of the respective S*-suffixes in T; hence, SA' determines 
the sorting of the S*-suffixes in T. Further, as at most every second suffix in T can be of type S*, 
the complete algorithm has worst-case running time T{n) = T(n/2) + 0{n) = 0{n), provided that 
the naming of the S*-substrings also takes linear time, which is what we explain next. 

The naming of the S*-substrings is similar to the inducing of the S-suffixes in the induced sorting 
algorithm (steps 2-4 above), with the difference that in step 2 we put the unsorted S*-suffixes 
into their corresponding buckets (hence they are only sorted according to their first character). 
Steps 3 and 4 work exactly as described above. At the end of step 4, we can assign names to the 
S*-substrings by comparing adjacent S*-suffixes naively until we find a mismatch or reach their 
end; this takes overall linear time. 
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3 Inducing the LCP- Array 



We now explain how the induced sorting algorithm (Sect. 12. 2p can be modified to also compute 
the LCP-array. The basic idea is that whenever we place two S- or L-suffixes Si-i and Sj-i at 
adjacent places k — 1 and k in the final suffix array (steps 3 and 4 in the algorithm), the length of 
their longest common prefix can be induced from the longest common prefix of the suffixes Si and 
Sj. As the latter suffixes are exactly those that caused the inducing of Si-i and Sj-i, we already 
know their LCP-value i (by the order in which we fill SA), and can hence set LCP[A;] to £ + 1. 

3.1 Basic Algorithm 

We now describe the algorithm in more detail. We augment the steps of the induced sorting 
algorithm as follows: 

1'. Compute the LCP-values of the S*-suffixes (see Sect. I3.3| l. 

2'. Whenever we place an S*-suffix into its S-bucket, we also store its LCP-value at the correspond- 
ing position in LCP. 

3'. Suppose that the inducing step just put suffix SsAfj]-! into its L-type c-bucket at position k. If 
'S'sA[i]-i is the first sufhx in its L-bucket, we set LCP[/!;] to 0. Otherwise, suppose further that 
in a previous iteration i' < i the inducing step placed suffix SsAfj/]-! at A; — 1 in the same c- 
bucket. Then if i' and i are in different buckets, the suffixes S'saIi] and S'sAfi'] start with different 
characters, and we set LCP[A;] to 1, as the suffixes SsA,[i]-i find SsAfj']-! share only a common 
character c at their beginnings. Otherwise {i' and i are in the same c'-bucket), the length i of 
the longest common prefix of the suffixes SgAfi] and SgAfi'] is given by the minimum value in 
LCP[i' + all of which are in the same c'-bucket and have therefore already been computed 
in previous iterations. We can hence set LCP[fe] to ^ + 1. 

4'. As in the previous step, suppose that the inducing step just put sufhx SsA\i]-i into its S-type 
c-bucket at position k. Suppose further that in a previous iteration i' > i the inducing step 
placed suffix S'sAfi'j-i at + 1 in the same c-bucket (if k is the last position in its S-bucket, 
we skip the following steps). Then if i' and i are in different buckets, their suffixes start with 
different characters, and we set LCP[/!; + 1] to 1, as the suffixes 5'sA[i]-i and S'sAii'j-i share only 
a common character c at their beginnings. Otherwise {i' and i are in the same c'-bucket), the 
length i of the longest common prefix of the suffixes SsA[i] and /SsAfj'] is given by the minimum 
value in LCP[i + l,i'], all of which are in the same c'-bucket and have therefore already been 
computed. We can hence set LCP[A; + 1] to £ + 1. 

3.2 Finding Minima 

To find the minimum value in LCP['i' + 1, i] or LCP[i + 1, i'] (steps 3' and 4' above), we have several 
alternatives. The simplest idea is to scan the whole interval from i' + 1 to i; this results in overall 
O(n^) running time. A better alternative would be to keep an array M of size a, such that the 
minimum is always given by M[c] if we induce an LCP-value in bucket c. To keep M up-to-date, 
after each step i we first set M[c] to LCP[i], and further update all other entries in M that are 
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larger than LCP[i] by LCP[i]; this approach has 0{na) running time. A further refinement of this 
technique stores the values in M in sorted order and uses binary search on M to find the minima, 
similar to the stack used by [7]. This results in overall O(nlgcr) running time. 

Yet, we can also update the minima in 0(1) amortized running time, as explained next. Let us 
first focus on the left-to-right scan (step 3'); we will comment on the differences to the right-to-left 
scan (step 4') at the end of this section. Recall that the queries lie within a single bucket (called c'), 
and every bucket is subdivided into an L- and an S-bucket. The idea is to also subdivide the query 
into an L- and an S-query, and return the minimum of the two. The S-queries are simple to handle: 
in step 3', only S*-suffixes will be scanned, and these are static. Hence, we can preprocess every 
S-type bucket with a static data structure for constant-time range minima, using overall linear 
space [11 Thm. 1]. The L-queries are more difficult, as elements keep being written to them during 
the scan. However, these updates occur in a very regular fashion, namely in a left-to-right manner. 
This makes the problem simpler: we maintain a Two-Dimensional Min-Heap [H Def. 2] A4c' for 
each bucket c', which is initially empty (no L-suffixes written so far). When a new L-suffix along 
with LCP-value £ + 1 is written into its c'-bucket, we climb up the rightmost path of A4c' until 
we find an element x whose corresponding array-entry is strictly smaller than i + 1 {Mc' has an 
artificial root holding LCP-value — oo which guarantees that such an element always exists). The 
new element is then added as x^s new rightmost leaf; an easy amortized argument shows that this 
results in overall linear time. Further, Aic' is stored along with a data structure for constant-time 
lowest common ancestor queries (LCAs) which supports dynamic leaf additions in 0(1) worst-case 
time [3]. Then the minimum in any range in the processed portion of the L-bucket can be found 



In the right-to-left scan (step 4'), the roles of the L- and S-buckets are reversed: the L-buckets 
are static and the S-buckets dynamic. For the former, we already have the range minimum data 
structures from the left-to-right scan (the 2d-Min-Heaps together with LCA). For the S-buckets, we 
now build an additional 2d-Min-Heap along with dynamic LCAs; this works because the S-buckets 
are filled in a strict right-to- left manner. 

What we have described in the preceding two paragraphs was actually more general than what 
we really needed: a solution to the semi-dynamic range minimum query problem with constant 
0(1) query- and amortized 0(1) insertion-time, with the restriction that new elements can only 
be appended at the end (or beginning, respectively) of the array. Our solution might also have 
interesting applications in other problems. In our setting, though, the problem is slightly more 
specific: the sizes of the arrays to be prepared for RMQs are known in advance (namely the sizes 
of the buckets); hence, we can use any of the (more practical) preprocessing-schemes for (static) 
RMQs in 0(1) worst-case time [HIS], and update the respective structures, which are essentially 
precomputed RMQs over suitably-sized blocks, whenever enough elements have arrived. 

^Note that it is important to use the Two-Dimensional Min-Heap rather than the usual Cartesian Tree for achieving 
overall linear time, for the following reason: Although the Cartesian Tree also has 0(1) amortized update-time for 
the operation "append at end;" it also needs to relink entire subtrees, rather than only inserting new leaves to the 
rightmost path 6 . For the relink-operation, no constant-time solutions exist for maintaining C'(l)-LCAs in the tree 
(not even in an amortized sense); the best solution we are aware of takes a{-,n) update time 18 , a(-, •) being the 
inverse Ackermann function. 
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3.3 Computing LCP- Values of S*-sufRxes 



This section describes how to compute the LCP-values of the suffixes in the sample set (step 
1' above). The recursive call to compute the suffix array SA' for the text T' (the text formed 
by the names of the S*-substrings) also yields the LCP-array LCP' for T' . The problem is that 
these LCP-values refer to characters Vi in the reduced alphabet [l,^'], which correspond to S*- 
substrings Ri in T. Hence, we need to "scale" every LCP-value in LCP' by the lengths of the actual 
S*-substrings that constitute this longest common prefix: a value LCP'[A;] refers to the substring 

^SA[fc] • ■ ■i'SA[fc]+LCP'[fe]-i of T', and actually implies an LCP-value of Yh^^^^'^ \Rsh{k\+i\ between 
the corresponding S*-suffixcs in T. 

A naive implementation of this calculation could again result in O(n^) running time, consider 
the text T = abab . . . ab. However, we can make use of the fact that the suffixes of T' appear 
lexicographically ordered in T': when "scaling" LCP'[/c], we know that the first m = min(LCP[A; — 
1], LCP [A;]) S*-substrings match, and can hence compute the actual LCP-value as 

LCP[fe]-l 

\Rs^.{k]+i 

i=0 



This way, by an amortized argument it is easy to see that each character in T contributes to at 
most 2 additions, resulting in an overall 0{n) running time. 

It is possible to stop the recursive LCP-calculation at a certain depth and use any other LCP- 
array construction algorithm on the remaining (sparse) set of sorted suffixes. 

3.4 Computing LCP-values at the L/S-Seam 

There is one subtlety in the above inducing algorithm we have withheld so far, namely that of 
computing the LCP-values between the last L-suffix and the first S-suffix in a given c-bucket (we 
call this position the L/S-seam). More precisely, when reaching an L/S-seam in step 3', we have 
to re-compute the LCP-value between the first S*-suffix in the c-bucket (if it exists) and the last 
L-suffix in the same c-bucket (the one that we just induced), in order to induce correct LCP-values 
when stepping through the S*-suffixes in subsequent iterations. Likewise, when placing the very 
first S-suffix in its c-bucket in step 4', we need to compute the LCP-value between this induced 
S-suffix and the largest L-suffix in the same c-bucket. (Note that step 4 might place an S-suffix 
before all S*-suffixes, so we cannot necessarily re-use the LCP-value computed at the L/S-seam in 
step 3'.) 

The following lemma shows that the LCP-computation at L/S-seams is particularly easy: 

Lemma 3. Let Si be an L-sujfix, Sj an S-suffix, and ti = c = tj (the suffixes are in the same 
c-bucket in SA). Further, let i > 1 denote the length of the longest common prefix of Si and Sj. 
Then 

Ti...i+t-i = = Tj jj^f^^i . 
Proof. Assume that t^+fc = c' = tj+fc for some 2 < k < i and c' 7^ c. Then if c' < c, both Si and Sj 



m—l m—1 

= X] l-^SA[fc]-hil + X] \^SA[k]+i 

t=n i=m 



already computed 
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are of type L, and otherwise (c' > c), they are both of type S. In any case, this is a contradiction 
to the assumption that Si is of type L, and Sj of type S. ■ 



In words, the above lemma states that the longest common prefix at the L/S-seam can only 
consist of equal characters. Therefore, a naive computation of the LCP-values at the L/S-seam 
is sufficient to achieve overall linear running time: every character ti contributes at most to the 
computation at the L/S-seam in the tj-bucket, and not in any other c-bucket for c ^ ti. 

4 Experimental Results 

We implemented the algorithm from the previous section in C and ran several tests on an AMD 
Athlon 64 processor, running at 2200 MHz with a 512KB L2-cache and 4GB of main memory. The 
basis of our implementation was Yuta Mori's linear-time C-implementation of the induced-sorting 
algorithm [18|, called sais-lite version 2.4.1 (http://sites.google.coni/site/yuta256/sais). We 
made the following implementation decisions: instead of calculating the LCP-values of the S*- 
suffixes recursively, we used a sparse variant of the ^-algorithm [10] immediately on the first level, 
which calculates the LCP-values of the S*-suffixes in overall linear time. For the inducing step, we 
used the simple 0(ncj)-variant described in Sect. 13.21 The resulting algorithm is called inducing 
henceforth. 

We compared our implementation to the following LCP-array construction algorithms: 

KLAAP: the original linear-time method for constructing LCP [12], implemented in a space-saving 
variant jl6j . 

$: the ^-algorithm of Karkkainen et al. [10], which is a clever variant of KLAAP that avoids 
cache-misses by reorganizing the computations. 

GO: the hybrid algorithm as described by [7J. It needs the Burrows- Wheeler- Transformation 
(BWT) for LCP-array construction, and computes small LCP-values naively, from which the 
larger LCP-values are deduced. 

G02: a semi-external variant of GO [7]. 

naive: for a sanity check, we also included the naive computation of the LCP-array (step through 
the suffix array and compare corresponding suffixes naively). 

We used the implementations from the succinct data structures library (sdsl 0.9.0) [7] wherever 
possible. All programs were compiled using the same compiler options (-ffast-math -09 -funroll- 
loops -DNDEBUG). 

We chose the test suite from 'http : //pizzachili . dcc.uchile . cl/" for evaluation, which is by 
now a de-facto standard. It includes texts from natural languages (English), biology (dna and 
proteins), and structured documents (dblp.xml and sources). Because the authors of [7] point out 
that the human chromosome 22 from Manzini's corpus (hs) is a particular hard case for some 
algorithms, it was also included. 
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The results are shown in Tbl. [TJ The first block of columns shows the running times for pure 
LCP-array construction. For KLAAP and ^, these times include construction of the inverse suffix- 
and the $-array, respectively, as they are needed for LCP-array computation. For GO and G02, 
the times for computing the BWT are not included; the reason is that in some cases the BWT is 
also needed for other purposes, so it might already be in memory. As inducing is inherently coupled 
with SA-construction flFJ, we could not measure its running times for pure LCP-array construction 
directly; the figures in column "inducing" of Tbl. [1] are hence obtained by first running the pure 
SA-construction (sais-lite), then the combined LCP- and SA-construction, and finally taking the 
difference of both running times. Measured this way, inducing takes always less time than all other 
methods. 

A fairer comparison of the algorithms is shown in the last three columns of Tbl. [U where the 
combined running times for SA- and LCP-array construction are given (for a selection of the 
best-performing LCP-algorithms). This is because all other methods for LCP-array construction 
are independent of the method for constructing SA, and can hence be combined with faster SA- 
construction algorithms. It is by now widely agreed that Yuta Mori's divsufsort in version 2.0.1 
is the fastest known such algorithm (http://code.google.eom/p/libdivsufsort/). Hence, for 
methods GO and naive we give the overall running times combined with divsufsort, whereas for 
inducing we give the overall running time of sais-lite, adapted to induce LCP-values as well. Further, 
for GO we also add the times to compute the BWT, as it is needed for LCP-array construction. 

Inspecting the results from Tbl. [H we see that inducing-l-sais-lite is usually the best possible 
combination, sometimes outperformed by naive-|-divsufsort. In fact, the naive algorithm is rather 
competitive (especially for small inputs up to 50MB), apart from the English text, which consists 
of long repetitions of the same texts (and hence has large average LCP). 

5 Conclusions and Outlook 

We showed how the LCP-array can be induced along with the suffix array. A rather ad-hoc im- 
plementation outperformed all state-of-the-art algorithms. We point out the following potentials 
for practical improvements: (1) As suffix- and LCP-values are always written to the same place, 
an interleaved storage of SA and LCP could result in fewer cache misses. (2) As the faster divsuf- 
sort is also based on induced sorting, incorporating our ideas into that algorithm could result in 
better overall performance. (3) Computing the LCP-values of the S*-suffixes recursively up to a 
certain (well-chosen) depth could be faster than just using the ^-algorithm on level 0, as in our 
implementation. 
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61.3 


130.3 


44.0 


proteins 


43.5 


38.6 


35.9 


46.7 


49.3 


16.2 


35.3 


48.7 


14.5 


85.7 


84.6 


64.9 


200MB 


dna 


104.4 


92.7 


46.1 


68.5 


51.0 


36.3 


75.9 


87.6 


32.7 


154.7 


126.9 


123.9 


english 


90.7 


80.9 


82.3 


104.3 


3190.5 


39.4 


68.9 


88.8 


31.6 


182.8 


3259.4 


128.2 


dblp.xml 


69.2 


64.6 


44.1 


64.3 


40.4 


31.1 


53.2 


63.6 


27.4 


124.7 


93.6 


94.7 


sources 


65.9 


62.0 


58.7 


79.9 


141.5 


29.5 


46.4 


65.3 


26.0 


131.1 


187.9 


94.8 


proteins 


91.6 


82.9 


82.0 


105.0 


124.2 


35.6 


76.5 


104.0 


30.8 


189.3 


200.7 


139.6 



As inducing is inherently coupled with SA-construction (sais-lite in our implementation), the running times for 



pure LCP-array construction were calculated by taking the difference of "inducing+sais-lite" and "sais-lite." 

Table 1: Running times (in seconds) for LCP- and suffix-array construction. The first block of 
columns shows the running times for pure LCP-array construction (for KLAAP and <1>, these times 
include construction of the inverse suffix- and the $-array, respectively). The second block shows 
the construction times of those arrays that need to be constructed before LCP: SA (always) and 
BWT (for GO and G02). The third block shows the overall running times for computing both SA 
and LCP for the best possible combinations of algorithms. 
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